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ABSTRACT 


Second-order and fourth-order two-dimensional advective equations 
developed by W. P. Crowley were examined to determine their applica- 
bility to atmospheric models. 

One second-order and two fourth-order forms were evaluated from 
their performances on a simple pattern advected. by a linear, divergent 
velocity field. The same equations were substituted for the advection 
term of a simple barotropic forecast model to determine their perform- 
ances on more general non-linear conditions. 

All forms of the equations remained stable in time and demonstrated 
the phase and amplitude characteristics predicted by Crowley. The 
fourth-order "advection" form gave best results. 

When substituted in the barotropic model, the fourth-order forms 
lead to improved trough and low-center movements, but RMSE was slightly 
larger than that resulting from second-order forms. The better RMSE 
of the lower order forms apparently resulted from their diffusive 


charateristics. 
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CHAPTER I 
INTRODUCTION 


The hydrodyamical equations governing fluid motion usually contain 
a non-linear advective or transport term which often plays an important 
role with respect to the temporal changes of the pertinent parameters. 
It is, therefore, desirable to approximate the advective term accu- 
rately when a numerical solution of the equations is sought. 

In order to investigate the accuracy of various finite-difference 
representations of the transport term, it will be assumed that some 
parameter, », is conserved in a two-dimensional flow field. This 


assumption may be expressed in Eulerian form as 


ad. V+ v9, (1.1) 


Q> 


— 
where V is the horizontal velocity and Vv is the two-dimensional del 
operator. 


When written in scalar form, eq. 1.1] becomes 


Boe. fu Sey B2, ae 


Q> 


where u and v are the x and y components, respectively, of the vector 


velocity. If the flow is non-divergent, eq. 1.2 may be expressed as 
o : 
oe = J(¥,o), (1.3) 


where ¥ is the stream function and J is the well-known Jacobian 
operator. 


An expression equivalent to eq. 1.2 15s 


ag. — |a(ou) , a(dov) gu, OV 


Although the analytical solution to eq. 1.4 is identical to that of 
eq. 1.2, the latter form leads to a different finite-difference expres- 
Sion possessing certain advantages which are discussed in Chapter II. 

For convenience, the advective equation as expressed by eq. 1.2 
will be called the "advection" form; equation 1.4 will be referred to 
as the "conservation" form. 

Although the analytical solution of a differential equation is 
exact, the finite-difference representation is only an approximation. 
When modern digital computers are utilized, the accuracy of the 
approximation is principally dependent on (1) truncation error, and 
(2) computational stability. Truncation error may be reduced (but 
never eliminated) by retaining high-order terms and by reducing the 
size of the grid-mesh. Elimination of instability error is more 
difficult because computational stability is dependent on the physical 
properties of the fluid as well as on the specific form of the finite- 
difference equation. 

Many finite-difference techniques have resulted from attempts to 
control instability while reducing truncation error. The purpose of 
this study was to evaluate by numerical experimentation the appli- 
cability to the atmosphere of high-order, two-dimensional advective 
techniques proposed by Crowley [1]. 

Examples of practical problems which require accurate advective 
techniques are (1) radioactive debris tracking, (2) air and water 
pollution studies, (3) estimation of movement and concentration of 
potentially lethal gases, and (4) numerical weather prediction. 
Numerical models which provide solutions to these problems often 


require a high degree of accuracy throughout long periods of time. 
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For example, radioactive debris projected into the stratosphere may 
travel around the earth several times before the combined effects of 
turbulent diffusion, velocity divergence and deformation, and fall-out 
reduce the concentration of the radioactive substance to an insignifi- 
cant level. 

A good advection scheme, therefore, must exhibit not only small 
errors in each time step, but also’ smal] cumulative errors after many 
steps. To determine the extent to which Crowley's equations satisfy 
this requirement, the initial experiment was designed to evaluate the 
performance of the techniques, given a known configuration of a 
conservative quantity and a known velocity field.’ From observations 
of the position, shape, and magnitude of the quantity after many hours 
of advection, statistics on effective phase speed, conservation of 
intensity and shape, and computational stability were obtained. 

Since diagnostic and prognostic models of the atmosphere must 
approximate real conditions in a satisfactory manner, the advective 
terms employed must demonstrate stability and accuracy under non- 
linear conditions. The characteristics of Crowley's equations, given 
non-linear conditions, were evaluated by comparing the results obtained 
from a simple barotropic forecast model in which appropriate forms of 


the various advective equations were utilized. 


CHAPTER II 
THE ADVECTIVE EQUATIONS 


The alternate forms of the advective equation will be developed in 
this chapter. Both result in a forward time step based upon centered 
space differences. Since the resulting finite-difference equations 
will be used eventually in the barotropic model, absolute vorticity, n, 
will be substituted for » in eas. 1.1 - 1.4. 

The following is an explanation of the notation used to indicate 
time and space differences: 

1. Superscript n: denotes the current hour. 


2. Subscript i: denotes the x-coordinate of the grid point, 
X = 1AX. 


3. Subscript j: denotes the y-coordinate of the grid point, 
y= aay. 


The development of equations 1.2 and 1.4 follows closely that of 


Crowley with minor differences in notation. 


THE ADVECTION FORM 


A Taylor series expansion in time of equation 1.2 in one dimension 





results in 
n+l n n Doee lt 
on (At)* 9°n 3 
N = 19) + At = i AyD Teo) ae be 
, UTM ae Te ae, (2.1) 
From equation 1.2, 
an. _ y an 
3 Ys (202) 


Differentiating the above equation with respect to time gives 


an. _ jau on, , 2 {an 
ote at ox at jax] |’ 
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If the order of differentiation is changed in the last term, and 


substitution is made from eq. 2.2, the result is 


an. {au an , a - u 3a) 











OX 

















377 = 2 on _ an{ au aah) dul 
at2 Ox? 9axtot Ox (28) 
Substitution of eq. 2.2 and 2.3 into 2.1 gives 
n 
naa = uae Sal" Cuged® én)" (as)* anl"fau _ y au 
r ‘ dX 2 dX¢ ; 2 Ox ee Ox ; (2.4) 





2 
If <2 and ao are evaluated to second order and the final term of 


eq. 2.4 is discarded, substitution of the finite-difference forms 





results in 
oe n f n n n 
LjuAt] 
n 72 yet in n | +4 21ax- ‘{n - 27 +7 . 
: " i+1 "es 1 41 0400 i+) (2n5) 


which is second-order accurate in space and will be referred to as the 
"second-order" form. Note, however, that if u is not constant, the 
equation is only first-order accurate in time. 

It is possible, although tedious, to develop, by similar methods, 
a form of eq. 2.5 which is fourth-order accurate in space. However, 
the "fourth-order" form so obtained remains only first-order in time. 


The finite-difference form is 
nt+1 n f oneoren fon n } n n 
LiUATi Ia} 
Ee 


¢ \¢ i en n n n 
- 54 (He I30n l6jn +7 | + In +7 } 
eee eee (G+) f-1 i+2 i-2 
cn n n n 1 
2in -n -In -7 | 
mee el te a 40) (2.6) 
Tale n (on n n n : 
ae (ta + in +7 
poles kee ge Woe? 4-2? 





THE CONSERVATION FORM 








The one-dimensional form of eq. 1.4 7 
an. , Bu _ aln 
at. SC «(ax ax (2573) 


The transport term, 3(ru)/a9x, can be considered to be a divergence of 
a flux. By Green's theorem, any change in the advected quantity with 
time is proportional to the differential of flux across the boundaries. 
If the boundaries are established at points half the distance between 


grid points, then the flux at one at oe boundaries is expressed by 


x ian 


yit 

1 
F = (mW) sya, * He {nl 
i+% it+% t : 
eG (2.8) 


i+ 


If the variation of n is assumed to be linear within the zone, then 


integration of eq. 2.8 produces 


~ 2 
uAt] uAt } 
Fay = hte |g Maen" tg + : 2a x] Poa nlf: (2.9) 


If the variation of n is described by a cubic equation centered on the 


point, i + 4%, the resulting equation is 


1 fuat) ° 
fut)” 
where the equation has been multiplied by 2 af, and uae 1s evaluated at 


| + ye 


When = is evaluated to second-order and the flux divergence is 


substituted tor the transport term, the second-order difference form 


Oi ed. 2./ iS 
nti noe a n n 
is At i At ‘ 
n =#=njlractu,.,- u, )ji- AF -F i), 


where the form of the fiux term is given by eq. 2.9, 
Evaluation of =< to fourth-order and substitution of flux diver- 
gence as computed by eq. 2.10, results in the fourth-order finite- 


difference conservation form 


n+l ne At ! 
a : Ws + Torx (8CUs- Yyig) - Uyyom U42))] 
-n n ) 
u ; ah (ae 
it% i-% 


STABILITY CHARACTERISTICS 
It is convenient to express the advective equations in terms of 
linear operators. Both forms of the advective equations described in 
this chapter have tne same general equation, 
n¢l f nj a] 


it. Teel gee 


A X é \ / 








where I represents the identity matrix and A represents the particular 
linear difference operator. Stability is, therefore, dependent on the 
eigenvalues of (i - A). Crowley [1] has shown that eqs. 2.5, 2.6, and 


Neo 
2.11 are stable for wee 1, and that eq. 2.12 is stable for ae < 


ax}! — 
ED. 


From Crowley's analysis of the amplitude and phase errors of the 
advection and conservation equations, the following generalizations 


are made: 


Ke 


Maximum amplitude damping and minimum phase retardation both 
oe 


occur at => 


The rate of decrease in damping and of increase in phase 
retardation is greater for values of woe larger than 0.7 than 
it is for values smalier than 0.7. 

Amplitude damping and phase retardation are serious only for 
wave lengths shorter than four grid intervals. 

The conservation form has slightly greater errors than the 
advection form in both amplitude and phase. 

The fourth-order forms are clearly superior to the lower order 


forms in both amplitude and phase. 
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CHAPTER III 
EXPERIMENTAL PROCEDURES 


Both forms of the advective equation were developed for one- 
dimensional systems. Two-dimensional systems may be constructed by 
employing Marchuk's method of fractional time steps [1,3]. The linear 


operator form of such a system is 


ol . ett 
n, , =/(1-A)(1-B)n ls .4 


1lsJ 


for the advection form, and 


n+1 


Nay -[(aky(a-an(a-ain |g (3.2) 


for the conservation form. An and Bn represent the fluxes in the x 
and y directions respectively, and k represents the velocity diver- 
gence term in the conservation equation. 

Using Marchuk's method, the final forms of the advective equations 


may be expressed as follows: 


* 8 
sr ga 
| 


Advection form: 


— 


bell - 

Mg? Coane (3.3) 
Conservation ea | a 

Nag 7 1K I-BI f 

ny - Ok )U-AIn (3.4) 


The symbol, — represents an intermediate value of the quantity at 
time n + 1, which results from the first one-dimensional step. Since 
each fractional step is, by itself, stable, the combination of the steps 


retains stability. 


2] 


EQUATION CHARACTERISTICS TEST 

The purpose of the initial experiment was to determine quantita- 
tively, the performance characteristics of the several equations under 
known conditions which have some similarity to the atmosphere. 

The grid system employed consisted of 3,969 points arranged in a 
Square matrix of dimensions 63 x 63. The interval, d, between grid- 
points was approximately two hundred nautical miles. In this grid, 
the point with coordinates (31,31) lies exactly at the center, so that 
distances from the center are proportional to [(31-i)2 + (31-42 lee 
For simplicity in this experiment, the surface was assumed to be 
flat; therefore no map factor was introduced. 

A velocity field was generated which was proportional to the 
distance from the center point and which contained no shear (solid 
rotation). A pure divergent component was then added to the rotation. 


The resulting equations are 


u=- My + Nx , 


and 


V Mx + Ny (325) 


where x = (i-31) and y = (j-31) are dimensionless coordinates. The 
first terms on the right are the rotational components and the second 
terms are the divergent components. M and N are arbitrary constants 
with units of meters per second, which determine the magnitude of the 
velocity components. 

The initial field of the quantity to be advected consisted of a 
sinusoidal perturbation centered on grid-point (31,19) and zeroes 


elsewhere. The equation of the perturbation is 


My? 10) ? sin(S)], R< 5, (ene) 


el A — 








where R = distance in grid intervals from point (31,19). Since d is 
approximately two hundred nautical miles, the perturbation has a 
maximum wave length of approximately two thousand nautical miles, 
which is of the same order of magnitude as the large scale vorticity 
disturbances. 

The computer program was designed so that the angular velocity and 
the divergence could be varied in sign and magnitude as well as the 
number and length of time steps. Values representative of the atmos- 
phere and forecast models were selected. The equation which relates 


the variables is 


2nR_d . MRS(At) (Se2) 
where 
ce = distance from grid center in grid intervals = 12, 
d = length of grid interval = 3.81 x 10° meters, 


S = number of time steps, 

At = length of time step, 

M = velocity factor (meters per second). 
Preliminary calculations determined that with one-hour time steps, 
approximately two hundred steps would be required to advect the 
quantity through 27 radians with a realistic wind. The necessary 


velocity factor is 


_ dnd » (6.28 )(3.81x10°) | 
M = Sat 2x0 3_6x103 = 3.32 m./sec. (3.8) 


This factor produces a linear velocity of approximately forty meters 
per second at the point (31,19). 
Initial experiments indicated that a suitable magnitude of N, the 


divergent factor, was 0.5 meters per second. 
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The scheme of the experiment was to advect the pattern through 1a 
radians with positive divergence, then to change the sign of the 
divergence for the remaining a radians. Assuming perfect advective 
properties, this method would return the pattern to its initial posi- 
tion and configuration after advection through 2n radians. Comparison 
of the initial and final patterns would show quantitatively errors in 
phase velocity and amplitude, as well as any computational instability 
which developed. 

For purposes of comparison, eq. 1.3 was also subjected to the 


above conditions. The finite-difference form of eq. 1.3 which was 





used was 
teat ge 5 on n n n 
Sy eee |) eer ey a 
a ne ee (3.9) 


which 18 second order in both space and time. 


THE BAROTROPIC EXPERIMENT 

The experiment just described subjects the advective equations to 
velocities which are linear in form and in comparison to the atmos- 
phere, rather idealized. In order to test the advective equations 
under more realistic conditions, a simple form of the barotropic 


forecast model was chosen as follows: 


Z 
vse 7 Eo 7 Pe (Gare 
where 
¥y = stream function 
y = constant = 4 
f = average value of coriolis parameter = 1.03126 x 107" sec.7! 
Z = average value of 500-MB height = 5574 m. 
Pe advective forcing function = V- wines 
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a 0 eee 


The form of the linear balance equation employed to compute the 
: A 


stream function was 


vey + a(vy.vF) = 2 vz, (3.11) 


However, the second term on the left in equation 3.11 was approximated 


by 


+ 


reducing eq. 3.1] to a Poisson-type equation. 

The scheme of the experiment consisted of equating the various 
forms of the advective equations to Fe This could not be done 
directly since F. = v- Vn, but by subtraction of In” from both sides 
of eq. 2.13, the desired quantity was obtained as follows: 


ntl n n 


n = laee=o- Agee= An, (3°12) 


The quantity An” is substituted for c in eq. 3.10. 

| Two different methods were used to integrate the stream function in 
time. When the Jacobian was used as the forcing function in the proa- 
ae cic equation, time integration of the stream function was accomplished 
by the leap-frog method. When the advert ton and conservation forms were 


used, time integration was done by forward steps, ij.e., 


n+l n : 
y = W + ay At : (3.13) 


which is equivalent to eqs. 3.1, 3-2, 3.3, or 3.4. 
apis 


PROGRAMMING PROCEDURES 


All numerical computation required by this study was accomplished 


A 


on the Control Data 6500 computer installed at Fleet Numerical Weather 


Central. FORTRAN, as modified for the CDC 6500, was the programmina 


ie x 


ZS 


language used. In addition to FORTRAN, several FNWC library sub- 
routines were utilized. 

No particular effort was expended to optimize the programs. 
Nevertheless, relative time requirements of the various methods should 
be significant. 

During the construction of the computer program to solve the 
finite-difference equations, several minor modifications of the 
advective equations were necessary. 

First of all, the boundaries of the grid necessitated the establish- 
ment of boundary conditions. The method chosen for the quantity, n, 
was simply to hold the boundary values constant. In the case of the 
fourth-order forms, the first interior value was computed using a 
second-order form. 

Since the flux value was evaluated at half arid-points, and the 
boundary values of n were not allowed to’ vary, the first computation 
of flux occurred between the boundary point and the first interior 
point. This was accomplished with a second-order form in the fourth- 
order equations. 

Similarly, the u and v components of the velocity field for the 
barotropic experiment were computed by second-order methods near the 
boundary, but in a different manner. 

Since the u and v components are derived from the gradient of ¥, 
jt was possible to use fourth-order methods along the boundaries which 
were perpendicular to the separate velocity components. For the 
boundaries which were parallel to the components, a second-order form 
was used for the first interior value. This particular method was 


necessary because the conservation form includes interpolated velocity 
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campenents and a diyergence term which requires the first derivatives 
of the velocity components. 

The remaining modifications were due to the fact that the flux 
divergence term was computed for the points 1 + % and j + %- Since 
the u and v components were computed for the interger points only, it 
was pecessary to interpolate tor the values at the intermediate points. 

This was accomplished by using a linear interpolation on the 
perpendicular boundaries and a second-order interpolation elsewhere. 


Fig. | shows the grid system employed in the conservation method. 


2] 


0 0 0 

0 i 0 a 0 r 0 
O 0 0 

0 : 0 + 0 AP O 
0 0 0 

0 z 0 5 0 a 0 
0 0 0 
a: +: aE 


Figure |. Schematic of grid used in the conservation 
method. Integer grid-points are represented 
by + 3; half-integer points by o. 
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CHAPTER IV 
EXPERIMENTAL RESULTS 


All contoured figures which appear in this chapter were produced 
by a California Computer Products, Inc , incremental plotter. The 
instructions to the piotter were computed by the CDC 6500 computer by 
a quadratic interpolation of grid-values Contour iabeis on figures 
applicable tu the jinear advection experiment are proportional to the 
actual vaiues The contour origin, 0, corresponds to an actual value 
of -1.0 and the contour interval corresponds to 1.5 units. 

To tacilitate the discussion of the various advective equations, 
the following symbols will be used: 


C4 denotes tou: th-order conservation 
A4 denotes fourth-order advection 

A2 denotes second-order advection 

Jl denotes common Jacobian 


EQUATION CHARACTERISTiCS TEST 

The experiment consisted of advection of a pattern through 2n 
radjans, about the center grid-point This waS accomplished in two 
hundred time steps of one hour each The initial pattern and the 
resulting patterns after advection through a and 2n radians were 
plotted Figs. 2-4 show the results when the fourth-order conserva- 
tion equation was used. Fig. 5 and fig 6 show the final results when 
the advection equations were used, 

Al} methods show good resuits, considering that the final pattern 
represents the cumulative errors in advection of the pattern once 


around the earth at a velocity representative of a moderate jet stream. 
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Figure 2. Contour plot of initial perturbation. Contour interval cor- 
responds to 1.5 units; contour label 06 corresponds to 0.5 units. 
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Figure 3. 





Plot of perturbation after 100 time steps by C4 method. 
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Figure 4, 





Plot of perturbation 


after 200 time steps by C4 method. 
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Figure 5. Plot of perturbation after 200 time steps by A4 method. 
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Figure 6. 


Plot of perturbation after 200 time steps by A2 
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method. 


There are, however, Significant quantitative differences. Table I 
shows the percentages of maximum error in phase and amplitude of each 
of the methods. 

For purposes of comparison, the experiment was repeated using the 
Jl equation. Prior to the completion of one hundred steps, instability 
near the grid boundaries compietely destroyed the pattern. This 
instability was eliminated by uSing one-half hour time steps. So 
that some reasonable comparison could be made between the methods, the 
number of time steps was doubled. The results are shown in fig. 7. 

A particularly interesting result of the latter test was that 
after advection over approximately equal distances using twice the 
number of time steps, the amplitude error of the Jl result was less 
than that of the A2 method. On the other hand, the phase error was 
larger than that resuiting from the A2 method. 

All of the results were qualitatively consistent with Crowley's 
theoretical computations which were summarized in Chapter III. The 
eccentricity of the resultant patterns was apparently due to the wave- 
length sensitivity mentioned by Crowley. Support for this statement 
was provided by a comparison of the initial and final positions of 
the outermost contour in the forward semicircles; they were nearly 
coincident in both fourth-order results. 

Although it is not apparent from the figures, small-amplitude 
noise was created by the advection and conservation methods. Nowhere 
did the amplitude of this noise exceed five percent of the magnitude 
of the perturbation. Except for the area close "behind" the pattern, 
the magnitude of the noise was less than one-half percent that of the 


perturbation. 
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PERCENTAGE OF MAXIMUM ERROR IN PHASE 


ADVECTIVE 
EQUATION A4 
PHASE 

ERROR 0.67% 
AMPLITUDE 

ERROR -25% 


TABLE I 


AND AMPLITUDE OF 
ADVECTIVE EQUATIONS 


C4 A2 
1.99% 3.45% 
- 34% -51% 
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- 48% 





teps 


Plot of perturbation after 400 one-half hour time s 
by Jl method. 
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Ftgure 7. 





A variation of the basic experiment was conducted to determine 
whether the length of the time step would have an effect. The time 
step was reduced to one-half hour and the number of time steps was 
doubled. The results obtained using the C4 method were nearly identical 
to those of the basic experiment. In this case, at least, any effect of 
reduction of the time step was apparently counteracted by the increased 
number of steps. 

In his conclusions, Crowley [1] surmised that given a divergent 
velocity field, the performance of the conservation equation would be 
Superior. He attributed the anticipated superiority to the explicit 
divergence term in the conservation form in contrast to the implicit 
divergence of the advection form. This experiment has not substan- 
tiated that supposition. 

Computation time is usually a consideration in the design of 
atmospheric models. Obviously, the high-order schemes treated in this 
study require more computer time than the mathematically more simple 
schemes. Due to the "purity" of the computation in this experiment, 
it is possible to estimate, with high confidence, the relative require- 
ments of the various forms by simply comparing the total computation 
time required. 

If the time required by the Jl method is considered to be unity, 


then the relative "cost" of the other methods is listed below. 


C4 = 5.14 
A4 = 3.74 
A2 = 1.98 


It is obvious that the cost of the high-order equations in a pure 


advective model is considerable, but when these equations are used in 
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an iterative model such as the barotropic forecast model, the percentage 


increase in total time is small. 


THE BAROTROPIC EXPERIMENT 

The 500-MB anaiyses used in this expertment were produced by the 
Fleet Numerical Weather Central operational upper-air analysis 
program Eight consecutive analyses beginning with 0000Z 15 January 
1968 and ending with 1200Z 18 January were processed. 

The choice of period was random except for the requirement that 
it be a winter period. As it turned out, it may have been an 
unfortunate choice because the upper air codes were changed on 
1 January 1968. A comparison of the average number of 500-MB reports 
successfully decoded during the period of study with the average number 
Successfully decoded during September 1968 results in the ratio, 
330/374. It is possible, therefore, that one or more of the analyses 
contained significant errors. 

Statistics on RMSE, pillow, variance of absolute vorticity, and 
mean absolute vorticity were computed for each map time. Contour maps 
of forecast values, and associated error maps, were produced for the 
following times only: 


¥200Z = 15 January 
O000Z 17 January 
1200Z 18 January 


The choice of these particular map times resulted from efforts to 
represent the extremes and mid-point of the period, and to avoid a 
multiple of twenty-four hours. 

In addition to the statistics based upon the results of the various 


advective equations, RMSE and pillow were computed from the operational 
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prognoses of FNWC. Contour plots of the FNWC operational forecasts 
corresponding to the three time periods listed above were also 
produced. 
Forecasts for twenty-four and forty-eight hours were produced. 
A one-hour time step was employed and no smoothing was applied to 
the forecasts produced by the advection and conservation equations. 
Early in the experiment, it was found to be necessary to smooth 
the stream fields when employing the Jl equation, due to short-wave 
instability near the boundary. The method adopted was to apply a 
heavy smoother to the first four interior rows of grid-points on the 
initial stream field and on the forecast stream field at the twenty- 
four hour point. The boundary points and other interior points were 


not smoothed. The equation of the smoother is 
Me + ¥;, + ¥ 
b= 1,j- 


+Y . 
my 62 Se at (4.1) 


This smoother was probably excessive, but no significant activity 
occurs south of 13°, the most northerly latitude of the smoothed points. 
The smoothing resulted in a significant reduction of O°» however. 

Table II shows the RMSE which resulted from the forecasts pro- 
duced by the C4, A4, A2, and Jl equations, and the FNWC operational 
forecasts. Means and standard deviations for the eight map times are 
also shown. 

Although RMSE is certainly not a’ conclusive measure of forecast 
skill, it is usually an indication of the relative accuracy of various 
methods. Table II suggested that in this experiment, skill was 
inversely proportional to the complexity of the advection equation 
employed in the model. This was indeed disturbing. Furthermore, it 


is surprising to note that all methods except C4, without smoothing, 
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Surpassed the skill of the FNWC operational model, which is considerably 
more complex than the model utilized in this experiment. 

To resolve this apparent paradox, the various prognostic charts 
and their associated error charts were compared with the verifying 
analyses. From this analysis the following conclusions were drawn: 


1. Closed lows existing on the analysis continued to deepen 
when C4 and A4 equations were used. When all methods 
continued the deepening process the C4 and A4 equations 
Showed the highest rate of deepening. Frequently, these 
lows showed filling on the verification maps. 


2. Movement of low centers was more rapid in the C4 and A4 
results, and usually more accurate when the low moved rapidly. 


3. Trough movements were significantly larger in the fourth-order 
results than they were in the results of the low-order equa- 
tions. Frequently, the trough movement forecast of the A4 and 
C4 equations exceeded the actual movement in the vicinity of 
the jet stream. 


4. None of the methods were successful in the development of a 
low-latitude cut-off low; a typical barotropic model failure. 


5. The southern extremity of all short-wave troughs, except those 
of FNWC, showed considerable lag. Considering this fact, the 
overall movement of troughs appeared to be the best in FNWC 
progs. Note, however, that FNWC has incorporated a modifica- 
tion tested by Lewit [4] to accelerate low-latitude troughs. 


From these conclusions, the inverse relationship of RMSE is at 
least partially explained. They do not, however, explain why RMSE of 
FNWC progs was large. Examination of the error charts showed that in 
general, FNWC error centers were smaller than those in the experimental 
methods. In the Pacific Ocean west of the Hawaiian Islands and over 
eastern Asia, however, FNWC errors were consistently larger than those 


in the experimental methods. A large positive error center in the 


42 





central Pacific and a moderate negative error over eastern Asia 
persisted on all FNWC charts. The Pacific error center appeared on 
charts of the experimental methods, but was of lower intensity. The 
negative error over Asia was absent in the experimenta! methods. No 
explanation was found for this inconsistency. 

Although generalizations were helpful in the explanation of the 
RMSE, they may not represent the characteristics of the model ina 
given forecast situation. A detailed analyses of the 15/1200Z map 
time was conducted to determine the qualitative and quantitative 
differences in performance of the models. 

The 15/1200Z map time was chosen because (1) the report count was 
the largest of all map times in the series, (2) the report count for 
the verification time, 16/1200Z, was well above average, and (3) the 
RMSE of the various methods closely approximated the mean values for 
all eight periods. 

Fig. 8 is a contour plot of the 15/1200Z analysis. It is typical 
of a rather low-index pattern with a nearly symmetrical blocking situa- 
tion in the central Pacific and deep troughs over the eastern United 
Brees. eastern Atlantic, and eastern Europe 

Fig. 9 is a plot of the 16/1200Z analysis, which provided the 
verification for the prognoses. It is apparent that little change 
has occurred in the Pacific hemisphere, but that the movement of 
Short waves in the Atlantic hemisphere has caused changes sufficient 
to exercise the forecast models. 

Figs. 10-14 are plots of the twenty-four hour progs produced by 
each of the methods. The most striking feature of these figures is 


the smoothness of the advection and conservation methods compared to 
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Figure 8, FNWC height analysis at 500-MB for 1200Z 15 January 1968. 
Contour interval is 60 meters. 
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Figure 9. FNWC height analysis at 500-MB for 1200Z 16 January 1968. 
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Figure 10. C4 24-hour 500-MB prog verifying 1200Z 16 January 1968. 
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Figure 11. A4 24-hour 500-MB prog verifying 1200Z 16 January 1968. 
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Figure 12. A2 24-hour 500-MB prog verifying 1200Z 16 January 1968. 
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Figure 13. Jl 24-hour 500-MB prog verifying 1200Z 16 January 1968. 
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Figure 14. 





FNWC 24-hour prog verifying 1200Z 16 January 1968. 
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the Jl method. This is a consistent feature. The smoothness of the 
FNWC prog is a result of repeated smoothing during the integration 
process. 
The trough over Labrador and associated ridges to the east and west 
appear to have been handled best by the C4 method. The FNWC prog is 
| nearly aS good in position and shape but has a much larger error in 
central pressure. This is a deepening situation. 

Large differences occurred in the handling of the low which, in the 
analysis, is located off the northeast Labrador coast. This low filled 
by verification time. The Jl method and FNWC excelled in this case. 
Note that the Crowley methods moved this feature northeastward and the 
A4 and A2 methods deepened it. Contrary to the generalization made 
earlier, the C4 method showed slight filling. 

The trough in the eastern Atlantic was handled rather well by all 
methods, but the J] method appears to have the best shape. Note that 
the A4, C4 and FNWC methods all move the trough too rapidly in the 
higher latitudes. 

The weak low located over southern Norway is of particular interest, 
since it is the only rapidly moving low on the chart. It can be fol- 
lowed to a position approximately 300 miles north of the Black Sea by 
verification time. The movement forecasted by the C4 method is 
clearly superior, but the central pressure is in error by more than 
120 meters, which is more than twice the error shown in the FNWC prog. 

The deep trough near the Caspian Sea was handled poorly by all 
methods. The fourth-order methods are particularly bad because of 
excessive deepening, especially the C4 method. Central Asia was 
handled similarly by all methods with only small differences in detail. 


The FNWC prog, however, appears to be overly smooth. 
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Although there appeared to be but small differences in the handling 

of the deep low north of Japan, the strong gradient produced quite 
different error patterns. The low moved northeast and deepened. 
The Jl, A4 and A2 methods moved it north. The FNWC prog resulted in 
movement slightly west of north. The C4 prog was clearly superior in 
the movement and deepening of this center, but also had an excessive 
northerly component. 

Since the western Pacific area was particularly difficult for the 
FNWC prog, it was interesting to compare the error charts of the FNWC 
prog with the C4 prog in this area. Fig. 15 and fig. 16 are computa- 
tions of prog minus verification, contoured at intervals of sixty 


meters. The differences are quite significant. 


All methods treated the blocking situation similarly. All failed 
to reduce the intensity of the low to the value shown on the verifica- 
tion. Note, however, that the low is positioned in a Sparse data area. 
It is quite possible that the center was more intense than indicated 
by the analysis. 

No significant differences in the handling of the Gulf of Alaska 
cyclone were apparent. 

The foregoing analysis of movements and intensities has indicated 
that no single method was always best. The C4 method usually gave best 
results when a system was intensifying or moving rapidly, but the FNWC 
prog, and often the low-order progs, appeared to give best results 
when the system was static or weakening. 

Investigation of the forecasts produced from the 17/0000Z and 
18/1200Z maps indicated that the above conclusions are applicable to 


other map times of this series. It does not follow, however, that 
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' Figure 15. C4 24-hour 500-MB prog error. Contours are prog minus 
| verifying analysis. Contour interval is 60 meters. 
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Figure 16. FNWC 24-hour 500-MB prog error. 
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these conclusions are applicable to all situations of the atmosphere. 
High-index situations, for example, might be troublesome, although the 
high-order methods appeared to excel when forecasting large movement. 
Weak, Summertime gradients might result in errors which are not 
apparent in this series. There is little doubt, however, that the 
advection and conservation equations remain stable in typical meteoro- 
logical wind fields. Moreover, the theoretical phase and amplitude 


characteristics appear to remain valid for these fields. 
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CHAPTER V 
SUMMARY 


Second-order and fourth-order advection equations have been 
developed for one dimension in conservation and advection forms. 
These equations have previously been shown to be stable for forward 


UAt < 1.0, and in the case of the fourth-order 


AX 
ee 
conservation form, when |~ | < 1.95. 
\ Ps 


Phase retardation was also shown to be a minimum and amplitude 





time steps when 


damping a maximum for van On 


somewhat better phase and amplitude characteristics than the 


The advection form demonstrated 


conservation form. 

Use of Marchuk's method of fractional time steps permitted the 
one-dimensional equations to be applied to a two-dimensional finite- 
difference scheme without loss of stability or modification of phase 
and amplitude characteristics. 

In this study a controlled experiment utilizing a divergent, but 
linear velocity field was conducted. A known pattern representative 
of a simplified atmosphere was advected through 2m radians by a 
velocity field of magnitude comparable to the atmosphere. Two hundred 
time steps were required. 

Positive, constant divergence was applied during the initial one 
hundred time steps; then negative divergence of equal magnitude was 
applied during the final one hundred steps. The object of the experi- 
ment was to determine how well the initial pattern was conserved. 
Crowley's theoretical stability analysis was substantiated by the 


experiment. 
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With a one-hour time step, computational instability of the 
Jacobian form destroyed the field before advection through w radians 
could be accomplished. Reduction of ust to one-half of the value 
used in the basic experiment and doubling the number of time steps 
eliminated the instability. It was found that the Jacobtan form had 
better amplitude characteristics than the second-order advection 
method, but that the phase characteristics were inferior. 

The several advection methods were substituted for the advective 
term in a simple non-divergent barotropic model. The purpose of this 
experiment was to test the characteristics of the methods in a realistic, 
non-linear situation. Operationaily produced 500-MB analyses were 
obtained from FNWC for initial data and for forecast verification. 

Twenty-four and forty-eight hour forecasts of the 500-MB height 
field were produced using each form of the advective equation on eight 
consecutive map times beginning at O000Z 15 January 1968. RMSE 
Statistics were compiled for each of the methods on each of the map 
times. Machine plotted contour maps of all forecasts were generated 
for the map times 15/i200Z, 17/0000Z, and 18/1200Z. Error charts were 
also plotted. 

Analysis of the RMSE statistics indicated that the relative skill 
of the methods was inversely proportional to the complexity of the 
method employed. Furthermore, all experimental methods except the 
fourth-order conservation form had jiower RMSE than did the FNWC 
forecasts. 

Subjective analysis of the prognoses and error charts indicated 
that the inverse relationship was principally due to a tendency of the 


high-order equations to deepen existing low pressure centers excessively. 


a7 


The poor performances of the FNWC progs was attributed to persistent 
errors in the western Pacific and eastern Asia. The cause of the FNWC 
errors is unknown. 

Detailed analysis of the 15/1200Z map time indicated that the 
fourth-order methods resulted in more rapid movement of troughs in 
strong gradient areas than occurred when lower order forms were used. 
Rapidly moving and deepening low centers were also handled better by 
the high-order methods, particularly the conservation method. Weakening 
low centers and static situations were generally handled better by the 
low-order methods. 

FNWC progs gave the best overall trough movement, primarily 
because of more rapid movement of the southerly portions. FNWC progs 
were at least equal to the lower order experimental methods when fore- 
casting for weakening cyclones. 

Crowley's equations were shown to be stable when subjected to the 
non-linear conditions of typical weather charts. The theoretical phase 
and amplitude characteristics appeared to remain valid under these 
conditions. 

No definitive conclusions concerning the relative skill of a 
barotropic model employing the high-order advective equations could 
be made. A large number of trials under widely varying initial condi- 
tions would be necessary. Such further testing appears to be warranted 
and desirable in view of the favorable results of the limited experi- 


ments conducted thus far. 
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